Energy Spectra of Quantum Turbulence: Large-scale Simulation and Modeling 
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In 2048 3 simulation of quantum turbulence within the Gross-Pitaevskii equation we demonstrate 
that the large scale motions have a classical Kolmogorov-1941 energy spectrum E{k) oc fc _5//3 , 
followed by an energy accumulation with E(k) ~ const at k about the reciprocal mean intervortex 
distance. This behavior was predicted by the L'vov-Nazarenko-Rudenko bottleneck model of gradual 
eddy-wave crossover [J. Low Temp. Phys., 153, 140-161 (2008)], further developed in the paper. 

PACS numbers: 25.dk,47.37.+q 



I. INTRODUCTION 

Hydr •adynamic turbulence (HD)i - loosely defined 
as a random behavior of fluids - remains the most impor- 
tant unsolved problem of classical physics, as was pointed 
out by Richard Feynman. 

Quantum turbulence (QT) - a trademark of turbu- 
lence in superfluid 3 He, 4 He and in Bose-Einshtein con- 
densates of cold atomic vapors 2 - - has added a new twist 
in the turbulence research shading light on old problems 
from a new angle. QT consists of a tangle of quantized 
vortex lines with a fixed core radius ao and a finite (quan- 
tized) velocity circulation k = h/M, where M is the 
proper atomic mass 2 -. The superfluid has zero viscosity, 
and in the zero-temperature limit, which is the simplest 
for theoreticians and reachable for experimentalists 3 -, the 
QT's Reynolds number, Re, is infinite. This brings (at 
least, the zero-temperature) QT to a desired prototype 
for better insight in the classical HD turbulence. 

The tangle of vortex lines in QT is characterized by 
a mean intervortex distance, I. For large i?-scale 
motions with i? 3> I the vortex tangles are better under- 
stood as bundles of nearly parallel vortex lines with mean 
curvature of about Ft£. For large scales the quantization 
of vortex lines can be neglected and QT can be considered 
as the classical one, in which the energy density in the k- 
space, E(k), is given by the celebrated Kolmogorov-1941 
(K41) law 4 .: 

E Kil (k) = Ce^k-^, (1) 
E(r) ee (\u(r)\ 2 )/2=jE(k)dk, 

confirmed experimentally and numerically-^. Here C ~ 1, 
e is the energy flux over scales, and E(r) is the energy 
density of turbulent velocity fluctuations per unit mass. 
Kelvin waves (KWs) are helix-like deformations of vor- 
tex lines with wavelength A: ao < A < I. Interactions of 
KWs on the same vortex line, but with different k ~ A -1 
lead to the turbulent energy transfer toward large k. This 
idea (Svistunov*) was developed and confirmed theoreti- 
cally and numerically by Vinen et al&, Kozik-Svistunov 



(KS)£ and L'vov-Nazarenko (LN)&. Two versions of KW 
spectrum where suggested in Refs. @Q : 

E KS (k) = C KS Ae 1/5 K 7 ^e- 8/5 k- J / 5 , KS, (2a) 
£ LN (fc) = C^Ae 1 / 3 /^- 2 / 3 ^ 4 / 3 ^- 5 / 3 , LN, (2b) 
A ee ln(£/a ) , 

Here C KS > 1, C LN ~ 1/tt (Ref. H), and * < 1 char- 
acterizes the ratio of I to the large-scale modulation of 
the vortex lines. Parameter A ~ 12 15 in typical ex- 
periments in 3 He and 4 He 3 -. The choice between Eqs. ^ 
is under intensive debates^—, which, however, has no 
principle effect on issues discussed in this paper. 

The nature of energy transfer and energy spec- 
trum is under intensive debates, too. Considering the 
incrtial (Re— > oo) energy transfer at the crossover scale 
k ~ L'vov-Nazarenko-Rudenko (LNR) pointed outi 3 - 
that for k ~ i^ 1 and A> 1 the KWs have much larger 
energy than the HD energy ([T]) at the same energy 
flux e. As the result LNR predicted a bottleneck en- 
ergy accumulation around k ~ £ . On contrary, KS 
suggested^ 4 - an alternative scenario due to possible dom- 
inance of vortex-reconnections in the energy transfer at 
k ~ without any energy stagnation. In Ref. [lfj] LNR 
predicted two thermal-equilibrium regions between the 
HD ([1]) and KW J2J energy-flux spectra: with equiparti- 
tion of the HD energy, E(k) oc fc 2 , followed by equiparti- 
tion of KW energy, E(k) ~const. 

The direct numerical simulations (DNS) of QT 
mostly use the Gross-Pitaevskii equation (GPE)i£,, which 
in dimensionlcss form is given by 



2idip/dt + V 2 il) = g\ip\ 2 ip 



(3) 



The macroscopic wave function ip(r,t) plays a role of the 
complex order parameter, and g is the coupling constant. 
The transformation ip = ^fp e %e maps Eq. ((3]) to the Eu- 
ler equation for ideal compressible fluid of density p and 
velocity u = V#, and an extra quantum pressure term. 

The numerical study of QT by GPE ((3]) has been re- 
ported in a few papers so far. Nore et al^ solved the 
GPE with resolutions up to 512 3 and observed that as the 
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quantized vortices became tangled, the incompressible ki- 
netic energy spectra seemed to obey the K41 law 
for a short period of time, but eventually deviated from 
it. Kobayashi and Tsubota^ solved the GPE on 256 3 
grid with an extra dissipation term at small scales, and 
showed the K41 law JT]) more clearly. Yepez et alS^ sim- 
ulated the GPE on grids up to 5760 3 by using a uni- 
tary quantum lattice gas algorithm. They also found a 
spectrum E(k) oc fc -5 / 3 and interpreted it as the K41 
law (TTJ of HD turbulence. However, due to the choice 
of initial conditions, their simulation should correspond 
to the pure KW region k > i~ x (thus supporting the 
LN-spectrum (J2bj) of KWs). 

In the present paper, we solved the GPE on the 
grids up to 2048 3 by parallelizing the simulation code on 
the Earth Simulator--. In contrast to Rcf. [2l| we fo- 
cused on HD- and crossover- regions, k < l~ x . 
First, we confirmed the K41-law (TT|) in the HD-region of 
about two decades long, which is wider than that of any 
previous work. 

Second, the visualization of vortices clearly shows the 
bundle-like structure, which has never been confirmed in 
GPE simulations on smaller grids. 

Third, we discovered a plateau in the crossover region, 
ki > 2tt, further explained as the KWs energy equipar- 
tition in the framework of the LNR's bottleneck models, 
which is revised here to account for the recently pre- 
dicted^ and numerically observed^ LN spectrum (|2b[) 
of KWs. 

We consider this correspondence as a support in favor 
of LNR bottleneck theory, understanding, nevertheless, 
that interpretation of numerical (or experimental) data 
with the help of a theoretical model on the edge of its 
applicability (here A ~ 1) is often problematic, being 
a question of experience, physical intuition and taste. 
Currently we cannot fully exclude the alternative KS- 
scenarioi^, even though it gives no energy stagnation for 
A ^> 1. More theoretical studies, numerical and labo- 



ratory experiments are required to fully understand the 
vortex dynamics in the crossover region of scales. 



II. NUMERICAL PROCEDURE AND RESULTS 

In DNS we follow techniques^ but extend the maxi- 
mum computational grid size from 256 3 up to 2048 3 . The 
initial state is prepared by distributing random numbers 
created inside a range from —Nira to Nira into the phase 
9(r) on selected points M 3 (M <C N) and interpolating 
them to make a smooth velocity field on all grid points. 
Here N is the total number of grid points and a is a 
control parameter for the initial energy input. Also, fol- 
lowing^ we add to the GPE an effective artificial en- 
ergy damping for small-scale motions by replacing in the 
Fourier transform of GPE i — > i+1 for k x , k y , k z > In/^, 
where £ ~ ao is the condensate coherence length. 

GPE conserves the total number of particles and the 
total energy (Hamiltonian) of the system^. We de- 
composi 
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the total energy density into the internal, 
E- m t = g(p~ l) 2 /4, the quantum, E qnt = \V^/p\ 2 /2, and 
the kinetic, E = p\u\ 2 /2, energy densities. The kinetic 
energy is decomposed into compressible and incompress- 
ible components, both of which are monitored. Two typ- 
ical spectra of the incompressible component are plotted 
in Fig. [TJ with corresponding vortex distributions. The 
plot illustrates a 512 3 run at times 3.8 and 7.8 in the 
left and right panels, respectively 3 ^. The time evolu- 
tion of the equation is calculated by a symplectic inte- 
gral method, and a typical pseudo-spectral method is em- 
ployed for the calculation of the kinetic energy term. The 
method is a standard one, which is known to guarantee 
sufficiently high accuracy for hydrodynamics simulations. 
In Fig. [TJ left, one finds that the major part of the energy 
spectrum fits the K41 law (fTJ) like ini2, but with the large 
incrtial interval. 

As expected, we also observed tangled vortex bundles 
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FIG. 2: Color online. Left: Simulation results of the incompressible energy spectra E(k£) normalized by fi 2 /(4m 2 £). Symbols: 
2048 3 (- • -), 1024 3 (-A-), 512 3 (-T-). A varies from A ~ 1.5 (slightly depending on time) for 512 3 to A ~ 2.2 for 2048 3 . 
Left: Dot-dashed (cyan) line: K41 "—5/3" scaling. A snapshot of vortex lines at the fully-developed turbulent state of 2048 3 
demonstrates the self-similarity of the bundle-vortex structure (see the dotted circles representing the zoom regions which 
vortex distributions are shown subsequently), typical for fully developed turbulence. 



clearly demonstrated in the insets of Fig. Q] showing a 
x-y 2D slice of the polarization field's color map, which 
is defined by summing up vortices (±1) inside plaque- 
ttes lying within a constant radius (= 32Aa?) from a grid 
point. On the other hand, Fig. Q] is a typical example of 
a considerably decayed state, in which the main feature 
is rather small vortex rings distributed almost equally 
inside the simulation cubic region. 

An important observation (Fig. [!} is a plateau-like re- 
gion for fc£ > 1.5 - a definite pile-up over the K41 spec- 
trum - a clear manifestation of the energy stagnation. 

The main numerical result of the present paper is 
Fig. [2] The left panel shows an inter-comparison of 
the incompressible kinetic energy spectrum E(k) among 
512 3 , 1024 3 and 2048 3 simulations. The K41 scaling ((TJ 
(shown as (cyan) dash-dotted lines) is extended to lower 
k range with the grid-size increase. This is the first clear 
demonstration of the classical K41 scaling characteris- 
tic for the normal-fluid turbulence but maintained in the 
large-scale range of the superfluid turbulence. The visi- 
ble extend of the K41 scaling on 2048 3 grid is much larger 
than that in all previous simulations. 

The right panel of Fig. [2] displays self-similar large 
structures of tangled vortices in the fully turbulent state: 
the large-scale vortex bundles in the maximum size, 
2048 3 , and smaller self-similar tangled structures inside 
this cubic region in the subsequent insets. 

Before discussion of these results we will revise shortly 
in the next Section the LNR model of the bottleneck 
crossover— to accont for recently predicted LN spectrum 
of Kelvin waves^. 



III. LNR MODEL OF THE BOTTLENECK 
CROSSOVER 

To find theoretically E(k) we, following LNR—, ap- 
proximate the superfluid motions as a mixture of "pure" 
HD- and KW-motions with the spectra E (k) = 
g{k£)E(k) and E KW (k) = [1 - g{kl)]E(k). Here g{ktj 
is the "blending" function, which was found in Ref. [la ] 
by calculation of energies of correlated and uncorrelated 
motions produced by a system of ^-spaced wavy vortex 
lines: 
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exp(a;) 



g(x) = 5o [0.321n(A + 7.5)x] , g (x) = [ . 

The total energy flux, also consisting of HD and 
KW contributions^, is modeled by dimensional rea- 
soning in the differential approximation. Hence, for 
k — > the energy flux is purely HD and thus Ek. oc 
k~ 2 \J E HD dE jdk. From the other side, for k — > oo the 
energy flux is purely KW and thus Ek oc [E ] 2 dE jdk. 
Important that in contrast to the Ref. [ill, where the 
physically irrelevant KS-spectrum of KWs (|2"a)) was used, 
we employ here the proper LN-spectrum (|2b[) that ac- 
counts for large-scale vortex-line modulations with short 
KWs£. The full equation for the total energy flux reads: 

2 



{\^k^g(k£)E{k) + l 



3 {^k 3 K l 2 [\ -g(k£)E(k)}}' 



(C LN A k) 



dk 



{E(k) 



dike) , l-g(k£) 
k 2 



+ 



kl 



Sk 



(4) 



Here £ HD (fc») 



E (fc.) , => kj ~ 6.64/ln(A + 7.5). In 
the inertial range the energy flux is constant e(fc) = e. 
Moreover, in the system of quantum filaments it is re- 
lated to the rms vorticity yj (\uj\ 2 ) ~ k/£ 2 via (l^l 2 ) = 
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A = 2 
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FIG. 3: Color online. Incompressible energy spectra plot- 
ted vs. kl and normalized by k 2 /£. The simulation results 
(the same symbols) and the LNR model for A = 10, 30, 100 
[dashed (black) curves] are brought together to the theoreti- 
cal [solid (black)] curve with A = 2 by superposing the K41 
(both simulations and model) and plateau regions (only for 
simulations). Dot-dashed (cyan) lines show different scaling 
asymptotics. 



2 / k 2 E(k)dk (see Refs. [10). This allows to find solu- 
tions of Eq. (j4} for different A as depicted in the Fig. [3] 
by (black) dashed and solid curves 3 ^. 

Four distinct scaling regions are evident (A> 1): 

I. k£ <C 1: E(k) and are dominated by the "pure" 
HD contributions, and the K41 law ([T]) is revealed. 

II. kl S> 1: E(k) and e k are dominated by the 
"pure" KW contributions, and one observes the LN spec- 
trum (|2bp of KWs with a constant energy flux. 

III. k < fc*: As explained above, for A ^> 1 the KW 
turbulence is much less efficient in the energy transfer 
over scales than its HD counterpart with the same energy, 
which leads to the (HD) energy accumulation with a level 
E(k) w E H °(k) > E K41 {k). For k < k*, both E{k) and e k 
are still dominated by HD contributions, but the energy 
flux is much smaller than the K41 estimate requires. This 
is like a flux-free HD system, thus, the thermodynamic 
equilibrium is expected with the equipartition of energy 
between the degrees of freedom: 3D-energy spectrum is 
constant, hence, the ID energy spectrum E (k) cx k 2 . 
This scaling is observed in Fig. 02 middle, for k < fc*. 
Think of a lake before a dumb, where the water veloc- 
ity being much smaller than that in the source river does 
not effect on the water level, which is practically horizon- 
tal. This interpretation of the energy bottleneck effect as 
"incomplete thermalization" (of only high k region) was 
suggested by Frisch et. al£^. 

IV. ^ Unexpectedly, we observe here almost 
a fc-independent ID-energy spectrum, E{k) ~ const, in- 
herent to the thermodynamic equilibrium of KWs. In 
the "pure" KW system, such a spectrum shows up for 
k ^> k*. However, in the region IV, the energy of the 



system is already dominated by the KW contributions, 
E(k) « E (k), while the energy flux is still dominated 
by the HD-motions 1 ^. Hence, this is almost a flux-free 
system of KWs, which is indeed found in the thermo- 
dynamic equilibrium with the ID energy equipartition: 

, ,KW 

E[k) = const. 

As one sees from Fig. [3] with the decrease of A the pile- 
up becomes less pronounced. For A = 2 the equilibrium 
HD region (III) almost disappeared, however the equilib- 
rium KW region (IV) is still well pronounced being much 
less sensitive to the value of A. 



IV. 



DISCUSSION AND SUMMARY 



A. Classical and quantum energy bottleneck effects 

The bottleneck effect in classical hydrodynamic tur- 
bulence is understood traditionall y 23 ^ 4 as a hump on a 
plot of compensated energy spectrum E(k)k 5 / 3 in the 
crossover region between inertial and viscous intervals. 
This is very general phenomenon reported in many nu- 
merical simulations and experiments of classical hydro- 
dynamic turbulence. For example, Yeung and Zhov*2£, 
Gotoh et al2£, Kaneda et al22 and Dobler et a\M found 
the bottleneck effect in their numerical simulations. Sad- 
doughi and Veeravalli 3 -^ studied the energy spectrum of 
atmospheric turbulence and reported the bottleneck ef- 
fect. Shen and Warhaft^, Pak et al^ 3 ., She and Jackson^ 
and other experimental groups also observed the bottle- 
neck effect in fluid turbulence. The bottleneck effect has 
been seen in other forms of turbulence as well, see e.g. 
Refs. [10-15] in Ref. H. 

To characterize the value of this effect one can intro- 
duce a "bottleneck magnitude" Mb n : the hump height, 
normalized by the plateau value of E(k)k 5 ^ 3 in the in- 
ertial interval. For example, in high resolution DNS of 
the classical hydrodynamic turbulence^ 3 ., shown in Fig.[fl 
left, its magnitude M bn ^ 0.34 for 512 3 DNS and de- 
creases with the resolution increase: Mb n — 0.31 for 
1024 3 and M hn ~ 0.25 for 2048 3 . Recent results^ based 
on the 4096 3 DNS confirm the statement that the bot- 
tleneck magnitude in classical turbulence systematically 
decreases with the DNS resolution increase (or equiv- 
alently, with the Taylor-Reynolds number ReA growth, 

and Mbn —> R- e A ° 4 as - ► °°)- 

Coming to comparison of the bottleneck effects in 
our modeling and numerical simulations we should no- 
tice that the LNR model accounts only for leading in 
A = ln(£/ao) term o 13 ' 15 . Moreover, it is based on the 
differential approximation for the energy flux, which is 
reasonable for vivid power-like behavior of the energy 
spectra, which exists only for A ^> 1, Fig. [3] Therefore, 
one expects that the LNR model is suitable for quan- 
titative analysis of experiments in 3 He and 4 He, where 
A ~ 12-^15, and can only qualitatively describe the sim- 
ulations presented here with A ~ 2. 

Nevertheless, the simulations clearly demonstrate in 



FIG. 4: Color online. Comparison of the compensated energy spectra, E(k) fc 5 ^ 3 (the same color-symbol code as in previous 
figures), and of the bottleneck magnitudes Mbn (dashed color arrows). Left: DNS data of classical hydrodynamic turbulence^ 
(cf. Figure 1 in Ref . [23h : modest values of M^ n i$ 1-34 decrease with the increase of the resolution. Right: Our DNS data of 
quantum superfluid turbulence: large values of Mhn increase with the resolution up to Mbn — 6.26. 



Fig. [3] the plateau that immediately follows the En- 
sealing which agrees with the LNR model predic- 
tion for A ~ 2 (Fig. [3]). The plateau broadens with the 
grid-size increase towards that of the LNR model curve 
(the earlier cutoff of the simulation data is due to the 
artificial dissipation) . The resolution of the current sim- 
ulations does not allow to resolve the KW-scaling (|2b[) 
with constant energy flux as it was done in^i, but the 
bottleneck is definitely there. 

To measure bottleneck magnitudes in quantum tur- 
bulence we re-plotted our data of Fig. [3J compensating 
E(k) by K41 prediction, i.e. multiplying by (fc£) 5 / 3 , see 
Fig. [4] right. One sees large humps with magnitudes Mb„ 
that increases with the resolution, reaching Mbn — 6.26 
for 2048 3 . Recall that in the classical turbulence Mb n 
is much (in about 20 times!) smaller and demonstrates 
opposite tendency with the resolution. 

We concluded that classical and quantum bottlenecks 
have completely different nature. Small magnitude of the 
bottleneck in classical turbulence is related to some non- 
locality of the energy transfer toward small scales that is 
slightly suppressed due to fast decrease of the turbulence 
energy in the dissipation range, while in quantum turbu- 
lence (at zero temperature) the essential bottleneck effect 
originates from the strong suppression of the energy flux 
in the Kelvin wave region. 

Indeed, Fig. [H right, demonstrates a good agreement 
between the QT DNS data and the LNR model predic- 
tion (that accounts for the flux suppression mentioned 
above) for A ~ 2, which improves with increasing the 
DNS resolution. The cutoff of the spectra for large k is 
a consequence of limited fc-space in the simulations. One 
predicts that with the further increase of the resolution 



the bottleneck magnitude can reach Mbn — 50 at A ~ 2 
and even much larger values for larger A. 

B. Summary 

In the paper we conclude that the observed essential 
bottleneck energy accumulation has definitely quantum 
nature (quantization of circulation) and can be com- 
pletely rationalized within the LNR model of gradual 
eddy- wave crossover, suggested in Ref. [Hj]. We con- 
sider this model as a Minimal Model of Quantum Tur- 
bulence that describes homogeneous isotropic turbulence 
in superfluids with energy pumped at scales much larger 
than the mean intervortex distance, and reveals reason- 
able (and even unexpectedly good) agreement with the 
simulations of Gross-Pitaevskii equation discussed here. 
The reason is that in the most questionable crossover 
region, the LNR model predicts a local thermodynamic 
equilibrium, where the energy spectra are universal and 
non-sensitive to the details of microscopic mechanisms of 
interactions, e.g. vortcx-reconnections, etc. 
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